Intermediate phase of the one dimensional half-filled Hubbard-Holstein model 



R.T. Clay and R.P. Hardikar 1 

department of Physics and Astronomy and ERC Center for Computational Sciences, 
Mississippi State University, Mississippi State MS 39762 
(Dated: February 2, 2008) 

We present a detailed numerical study of the Hubbard-Holstein model in one dimension at half 
filling, including full finite-frequency quantum phonons. At half filling, the effects of the electron- 
phonon and electron-electron interactions compete, with the Holstein phonon coupling acting as an 
effective negative Hubbard onsite interaction U that promotes on-site electron pairs and a Peierls 
charge-density wave state. Most previous work on this model has assumed that only Peierls or 
U > Mott insulator phases are possible at half filling. However, there has been speculation that a 
third metallic phase exists between the Peierls and Mott phases. We present results confirming the 
intermediate metallic phase, and show that the Luttinger liquid correlation exponent K p > f in this 
region, indicating dominant superconducting pair correlations. We explore the full phase diagram 
as a function of onsite Hubbard U, phonon coupling constant, and phonon frequency. 
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Electron-phonon (e-ph) interactions can give rise to a 
number of interesting effects in low-dimensional materi- 
als, including superconductivity as well as charge-density 
wave and insulating phenomena. Frequently these ma- 
terials feature strong electron-electron (e-e) interactions 
as well, leading to very rich phase diagrams that com- 
bine lattice, charge, and spin (magnetic) orderings. We 
will focus specifically on materials where the electrons 
are coupled to localized vibrational modes, which may 
be of relatively high frequency. This type of e-ph in- 
teraction is most studied in molecular crystal materials, 
including the quasi-one- and quasi-two-dimensional or- 
ganic superconductors Q and fullerene superconductors 
[2j. In all of these materials, a fundamental question is 
whether the effects of e-e and e-ph interactions compete 
or cooperate with each other. In this Letter, we examine 
this issue within one of the most basic models. We find 
that despite the two interactions each separately favor- 
ing insulating states, together they can mediate an un- 
expected metallic phase with superconducting (SC) pair 
correlations. 

The model we consider is the one-dimensional (ID) 
Hubbard-Holstein model (HHM), with Hamiltonian 

H = -t^{c] +la c 3 ^ + h.c) + uY^ n 3A n 3,l (!) 
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where Cj a (cj )Cr ) are fermionic creation (annihilation) op- 
erators for electrons on site j with spin a, bj (bj) are 
bosonic creation (annihilation) operators for phonons at 
site j, and n^ a — cj a Cj :<T . The dispersionless phonons 
have frequency u> and are coupled to the local electron 
density with coupling strength g 3]. U is the Hubbard 
on-site e-e interaction energy. All energies will be given 
below in units of the hopping integral t. 



The properties of the ID 1/2-filled HHM are well un- 
derstood in two limits: the static or lu — > and the 
u> — > oo limit 0. First, in the static limit, the ground 
state is Peierls distorted for any nonzero e-ph coupling 
g = + . As the phonons couple to the electron density, 
the Peierls state is a 2kp charge-density wave (CDW) 
consisting of alternating large and small site charges. In 
the u> — > oo limit, the retarded interaction between elec- 
trons mediated by the phonons becomes instantaneous 
in imaginary time, and the phonons may be integrated 
out. This leads to an effective renormalized Hubbard in- 
teraction U C B = U — 2g 2 /lu. While strictly at lo — > oo the 
Peierls state cannot occur, for finite lu the Peierls state 
may again occur, although the mapping to an effective 
negative U is no longer exact. However, based on the 
lo — > oo mapping, it has been believed that the ground 
state phase of the HHM may be determined through U e g. 
U e ff should correspond to the U > Hubbard model, 
which in ID at 1/2 filling has a finite charge gap for any 
U > and no spin gap. We shall refer to this state as the 
Mott state. On the other hand, if U c s < 0, e-ph inter- 
actions dominate over e-e interactions, and the ground 
state is Peierls CDW distorted. 

Hence at ui = 0, the Peierls distortion occurs uncondi- 
tionally, while at lu = oo there is no Peierls state. Much 
less is known in the intermediate u) region. In the U = 
model, the Peierls state may be viewed as a traditional 
band insulator for small w, and as bipolaronic insula- 
tor composed of tightly bound pairs in the large lu limit, 
with a crossover between these pictures for lo ~ t 0|. 
If lo > 0, the Peierls distortion only occurs for e-ph cou- 
pling g above some critical value g c F° r 9 < 9c the 
ground state is then metallic at 1/2 filling (at U = 0). For 
[7 = 0, g^has been calculated using a functional integral 
method jfj and also via the density matrix renormaliza- 
tion group (DMRG) [2j. For intermediate U as well as co, 
far less is known. However, it has been recently proposed 
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FIG. 1: U = 0, U! = 1 results for long-wavelength charge (open 
symbols) and spin (filled symbols) structure factors versus g 
for periodic systems of N = 16 (diamonds), N = 32 (circles), 
and N = 64 (squares) sites. Statistical errors are smaller than 
the symbols. The inset shows finite-size scaling of the critical 
coupling g c (indicated by arrow) where K p — 1. 

that a metallic ground state exists intermediate between 
the Mott and Peierls states at 1/2 filling, i.e. when U e g is 
close to zero 8] . This metallic phase occurs for intermedi- 
ate uj and hence cannot be predicted from the small and 
large u limits. Our goal in this Letter is to confirm this 
metallic state and investigate its properties. Its existence 
is perhaps not surprising given the known existence of a 
metallic phase for g < g c at U — 0; this region of the 
phase diagram continues to exist for finite U. We further 
show that this metallic region persists for a substantial 
range of parameters provided u> is not too small. 

Hamiltonian Eq. ^ is difficult to analyze in the inter- 
mediate coupling region due to the presence of both elec- 
trons and phonons. The numerical method we use is the 
Stochastic Series Expansion (SSE) quantum Monte Carlo 
method with directed loops |9j . SSE is a powerful method 
for non-frustrated quantum spin systems or ID electron 
lattice models where no sign problem occurs. Impor- 
tantly, there are no approximations in the method besides 
finite system size and temperature. Electron-phonon in- 
teractions have been incorporated in SSE for both spin 
models ^3 an d electron models 01- As m these refer- 
ences, we treat the phonons in the occupation-number 
basis. The number of phonons per lattice site is un- 
bounded in the thermodynamic limit, but for a finite 
system at a finite temperature, the number of phonons 
may be truncated. We choose this truncation in a similar 
manner as the truncation of sequence length in the SSE 
method 9J: in the equilibration phase of the calculation, 
the phonon truncation is increased to exceed the cur- 
rent number of phonons on any given lattice site when- 
ever necessary All results shown below used periodic 
lattices of ./V sites, with inverse temperatures of at least 
f3/t = 2N and phonon cutoffs of up to 30 phonons per lat- 



tice site. In this Letter we focus on results for the HHM 
model, and details on the SSE implementation will be 
published separately. Our code was checked extensively 
against Lanczos exact diagonalization results for several 
different observables. We also implemented the quantum 
parallel tempering algorithm |l2l | , where different proces- 
sors of a parallel computer are assigned different model 
parameters (U, g, and lu). A Metropolis probability is 
then computed to switch configurations between proces- 
sors with adjacent parameters. As in reference Il2t we 
found this technique essential in obtaining smooth data 
across quantum phase transition boundaries. 

The low energy properties of any ID gapless interact- 
ing electron model may be mapped to an effective con- 
tinuum model, or Luttinger Liquid (LL). The properties 
of the LL, and in particular the decay with distance of 
different correlation functions are then described by two 
correlation exponents, K p for charge properties, and K a 
for spin 0| . K p values greater than 1 indicate domi- 
nant attractive SC correlations (no SC long-range order 
is possible in strictly ID systems), while K p < 1 corre- 
sponds to repulsive charge correlations. For models with 
spin-rotation symmetry, K a is always equal to 1. K a = 
then indicates the presence of a spin gap. These expo- 
nents are most easily computed via the SSE data from 
the static structure factors: 

S »Aq) = ^£e i,0_fc) <KT ±^l)Kt ±n H )) (2) 

K p and K a are then proportional to the slope of the cor- 
responding structure factor in the long wavelength limit 
g^O Hi: 

Kp, a = — s p , a (q -» 0) (3) 

nq 

These values once finite-size scaled may then be used to 
determine the quantum phase boundaries. A second ob- 
servable we will use are the charge and spin stiffnesses, 
p c and measured in the SSE method via the winding 
number [12j. A zero stiffness indicates a gap in the cor- 
responding sector, while nonzero p indicates no gap. We 
also verified directly that charge-charge (spin-spin) cor- 
relation functions showed staggered order in the Peierls 
(Mott) phases. 

We first present results for U = and uj = 1. 
Fig. ^ shows the slope of charge and spin structure fac- 
tors irSp, a (qi)/qi evaluated at the smallest wavevector 
qi = 2n/N, plotted versus e-ph coupling g. For g — > 
both Kp and K a tend to exactly 1, as required for free 
electrons. Finite size effects become very small in this 
limit, consistent with the expected vanishing of logarith- 
mic corrections when K p = K a = 1 0, ^| . At a critical 
coupling g c , K p crosses 1 and tends to zero, marking the 
transition to the Peierls state. From finite-size scaling 
of 16, 32, and 64 site systems vs. 1/JV (see Fig. ^ in- 
set), we find the critical coupling g c = 0.66 ± 0.01, with 
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FIG. 3: Charge (a) and spin (b) stiffnesses versus g for (7=2, 
uj=5. Open (closed) symbols are for A r =32 (16) site lattices. 
Arrows mark the transition points g c \ and g C 2 determined 
from Fig.Hb). 



FIG. 2: Long-wavelength charge (open symbols) and spin 
(filled symbols) structure factors versus g. (a) U = 2, uj =0.5; 
arrow indicates finite-size scaled g c (b) (7=2, u=5; arrows 
indicate finite-size scaled transitions g c \ and g a 2- 



the uncertainty estimated from the linear fit. A previous 
DMRG study, which did not attempt finite-size scaling 
found g c « 0.8 |7(. For g < g c , the LL exponent K p > 1, 
and K a scales to zero with increasing N. Therefore, for 
U = and g < g Cl the ground state is metallic, with 
dominant SC pair correlations and a spin gap. 

For any finite U at g = and 1/2 filling, the domi- 
nant ground state correlations of Eq. ^ arc spin-density 
wave (SDW). In this ground state K a — 1, and K p = 
indicating an insulating state with no spin gap. As this 
Mott state is not equivalent to the metallic state found 
for g < g c in Fig. ^ we may then expect the possibility of 
three different phases for the HHM with U > and g > 0: 
Mott, Peierls, and metallic SC. Fig. [21 shows the charge 
and spin structure factors again versus g with [7 = 2. We 
first discuss results for small phonon frequency oj = 0.5 
shown in Fig. Ufa). We find that for small g, K a tends 
toward 1 as system size increases. At a critical coupling 
g c , K a becomes less than 1, indicating the opening of 
a spin gap and the transition to the Peierls state. The 
charge exponent K p scales towards zero on both sides of 
the transition, developing a peak at g — g c . Similar re- 
sults are found for the 1/2-filled ID extended Hubbard 
model (EHM), where K p also peaks at the transition be- 
tween CDW and bond-order wave (BOW) phases [l|. 



The finite value of K p on the boundary indicates that the 



transition is of continuous nature 



l ary : 
.16]. 



Like the CDW- 



(BOW)-SDW transition in the 1 /2-filled EHM, we find 
similar behavior for small uj in the HHM as U increases: 
the value of K p at the transition decreases with (7, indi- 
cating that for large U the transition is first order rather 
than continuous. The correspondence is not surprising, 
as the effect of a nearest-neighbor interaction V ^> U can 
be viewed as an effective negative ?/[l4j- 

Charge and spin response typical for large uj are shown 
in Fig. |2Ib), here shown for (7 = 2 and uj = 5. Again, 
for small <?, we find K a tends toward 1 as system size 
increases, and K p tends toward 0, consistent with the 
Mott state. However, at a critical coupling g c \ K p ex- 
ceeds 1 and K a becomes less than one, indicating the 
dominant-SC state found for (7 = 0, g < g c . The value 
of g c \ is very close to the expected value where (7 er r = 0. 
Increasing g further, for g = g C 2 K p becomes less than 
1, indicating the opening of a charge gap and the Peierls 
state. In the intermediate region g c \ < g < g C 2, the 
properties of the model are identical to the (7 = model 
for g < g c as seen in Fig. ^ i.e. metallic with dominant 
SC correlations. While the second transition point g c2 is 
finite-size dependent, finite-size effects are very weak at 
the first transition since K p = K a = 1. It is then clear 
that due to the crossing of K p and K a curves at exactly 
1 at g = g c \ , a region of K p > 1 must exist for g > g c \ . 

To confirm the structure factor results indicating two 
transitions at g c \ and <? C 2, we present the charge (p c ) 
and spin (p s ) stiffnesses in Fig. As with the structure 
factors, the presence of two transitions may only be de- 
tected from p c and p s after finite-size scaling has been 
performed. For clarity, we plot only two system sizes in 
Fig. El In Fig. 01a) we see that for small g, p c decreases 
with system size, indicating a charge gap. p s is nearly 
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FIG. 4: Evolution of the intermediate phase versus u. (a) 
w=0.5 (b) cj=1 (c) u>—5. Open circles indicate the boundary 
p c i between Mott and SC phases, filled diamonds the bound- 
ary gr C 2 between SC and Peierls phases, and open squares the 
boundary between Mott and Peierls phases. Dashed line in- 
dicates U — 2g 2 /uj. 

constant with system size, indicating no spin gap. In the 
intermediate phase, this is reversed, with p c remaining 
constant or increasing with system size, and p s decreas- 
ing with system size. This again indicates a spin gap but 
no charge gap in the SC region. Finally, for g > g C 2, both 
stiffnesses go to zero in the Peierls state. 

Fig.0]shows the evolution of the intermediate phase as 
a function of to, with phase boundaries determined from 
finite-size scaling of 16, 24, and 32 site structure factor 
data. In all cases, the SC phase exists for g < g c exactly 
at U — 0. As U is increased at fixed uj, the SC region 
then shrinks. For small enough cj, we see two different 
sequences of phases, either Mott-SC-Peierls for small U 
(as in Fig. Efb)), or Mott-Peierls (as in Fig. for 
large U. Numerically it is difficult to precisely deter- 
mine point where all three phases meet, but we find that 
the SC phase disappears at U « 1 for uj = 0.5, and at 
U k 2 for uj = 1.0. In Fig. 0fa) and (b), we have not 
plotted points for the Mott/Pcierls boundary for large 
U , as in this region the transition becomes strongly first 
order, making exact determination of the boundary dif- 
ficult. However, the transition appears to remain close 
to g 2 /uj = U/2. For uj = 5, we were not able to access 
large enough U to determine the upper cutoff U neces- 



sary to suppress the SC state, but the intermediate phase 
appears to persist up to at least U « 7 for uj — 5. In all 
cases, we find the line U c s = to very accurately predict 
either the Mott-SC or Mott-Peierls boundary, but not 
the SC-Peierls boundary. Our phase diagram is slightly 
different from reference |3, where the metallic region was 
found to extend to either side of the U c s = line. 

In conclusion, we have shown that a metallic SC region 
exists intermediate between Peierls and Mott phases in 
the ID 1/2-filled HHM. While in the 1/2-filled model 
with small phonon frequency uj < t the intermediate 
phase occupies a relatively small region of the phase di- 
agram, we expect that in the presence of doping (in fact 
most of the organic SC's are 1/4- filled 0), the size of 
this region will be greatly enhanced. Furthermore, while 
the SC pairing found here exactly at 1/2- filling consists 
of on-site pairs or bipolarons, with doping both onsite 
and nearest-neighbor pairing may be mediated by the 
Holstein phonons|l7j. We are currently exploring these 
possibilities in the doped model. 
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